library(sf)
## Warning: package 'sf' was built under R version 4.1.3
## Linking to GEOS 3.9.1, GDAL 3.2.1, PROJ 7.2.1; sf_use_s2() is TRUE
library(sp)
## Warning: package 'sp' was built under R version 4.1.3
library(gstat)
## Warning: package 'gstat' was built under R version 4.1.3
library(raster)
library(tidyverse)
## Warning: package 'tidyverse' was built under R version 4.1.3
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v ggplot2 3.3.6 v purrr 0.3.4
## v tibble 3.1.6 v dplyr 1.0.8
## v tidyr 1.2.0 v stringr 1.4.0
## v readr 2.1.2 v forcats 0.5.1
## Warning: package 'ggplot2' was built under R version 4.1.3
## Warning: package 'tidyr' was built under R version 4.1.3
## Warning: package 'readr' was built under R version 4.1.3
## Warning: package 'dplyr' was built under R version 4.1.3
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x tidyr::extract() masks raster::extract()
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
## x dplyr::select() masks raster::select()
library(spacetime)
library(stringr)
library(readr)
library(leaflet)
library(gridExtra)
##
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
##
## combine
library(gifski)
## Warning: package 'gifski' was built under R version 4.1.3
# to make the spherical geometry errors go away
sf::sf_use_s2(FALSE)
## Spherical geometry (s2) switched off
# AIR QUALITY
## corrected air quality
AQ_df = read.csv("../intermediary_outputs/corrected_AQ_data.csv")
## sensor info
sensors = read.csv("../intermediary_outputs/sensor_data_full.csv")
# BOUNDARIES
# Boulder boundaries
BO_CO = st_read("../GIS_inputs_destruction_fireboundary/Boulder_county_munis/Municipalities.shp")
## Reading layer `Municipalities' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\Boulder_county_munis\Municipalities.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 10 features and 10 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -105.525 ymin: 39.91405 xmax: -105.0528 ymax: 40.23771
## Geodetic CRS: WGS 84
(fire_counties = BO_CO %>% dplyr::filter(ZONEDESC == "Louisville" | ZONEDESC == "Superior" | ZONEDESC == "Broomfield" | ZONEDESC == "Lafayette" | ZONEDESC == "Boulder"))
## Simple feature collection with 4 features and 10 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -105.3014 ymin: 39.91405 xmax: -105.0528 ymax: 40.09448
## Geodetic CRS: WGS 84
## OBJECTID ZONECLASS ZONEDESC LASTUPDATE LASTEDITOR
## 1 458263 XLO Louisville 2016-03-02 WTODACHEENE
## 2 458380 XBO Boulder 2021-10-06 WTODACHEENE
## 3 458382 XLA Lafayette 2019-01-18 WTODACHEENE
## 4 458443 XSU Superior 2021-03-23 WTODACHEENE
## REG_PDF_UR Shape_STAr Shape_STLe ShapeSTAre ShapeSTLen
## 1 http://www.louisvilleco.gov 0 0 225645149 155909.46
## 2 http://www.bouldercolorado.gov/ 0 0 775961674 448684.76
## 3 http://www.cityoflafayette.com/ 0 0 263931665 239628.72
## 4 http://www.superiorcolorado.gov 0 0 105234151 74391.85
## geometry
## 1 MULTIPOLYGON (((-105.1824 3...
## 2 MULTIPOLYGON (((-105.2113 4...
## 3 MULTIPOLYGON (((-105.1063 4...
## 4 MULTIPOLYGON (((-105.1779 3...
boulder_precincts = st_read("../GIS_inputs_destruction_fireboundary/Unincorporated_Boulder/Unincorporated_Boulder.shp")
## Reading layer `Unincorporated_Boulder' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\Unincorporated_Boulder\Unincorporated_Boulder.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 129 features and 8 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -105.4297 ymin: 39.91371 xmax: -105.0528 ymax: 40.21256
## Geodetic CRS: WGS 84
broomfield_precincts = st_read("../GIS_inputs_destruction_fireboundary/Broomfield_Precincts/Precincts.shp")
## Reading layer `Precincts' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\Broomfield_Precincts\Precincts.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 61 features and 9 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 3093759 ymin: 1202697 xmax: 3150881 ymax: 1259386
## Projected CRS: NAD83(HARN) / Colorado North (ftUS)
westminster_city = st_read("../GIS_inputs_destruction_fireboundary/Westminster_CityLimits/CityLimits.shp")
## Reading layer `CityLimits' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\Westminster_CityLimits\CityLimits.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 34 features and 18 fields
## Geometry type: POLYGON
## Dimension: XY
## Bounding box: xmin: -105.1654 ymin: 39.81811 xmax: -104.9877 ymax: 39.96854
## Geodetic CRS: WGS 84
# make all the precincts have the same columns (we don't need a lot of the data)
boulder_precincts <- boulder_precincts %>%
select(c(OBJECTID, geometry, SHAPEarea, SHAPElen))
broomfield_precincts <- broomfield_precincts %>%
mutate(SHAPEarea = Shape__Are,
SHAPElen = Shape__Len,
OBJECTID = OBJECTID + 1000) %>% # add 1,000 to the broomfield precincts to get rid of overlapping labels
select(-c("GlobalID", "GIS_ID", "PRECINCT_N", "MAP_COLOR", "USER_COMME", "LAST_UPDAT", "Shape__Are", "Shape__Len"))
westminster_city <- westminster_city %>%
mutate(SHAPEarea = ShapeSTAre,
SHAPElen = ShapeSTLen) %>%
select(c(OBJECTID, geometry, SHAPEarea, SHAPElen))
# save proj string for fire-affected counties to use for other transformations
prg = raster::crs(fire_counties,asText=TRUE)
# all boundaries joined together
all_bounds <- boulder_precincts %>%
rbind(st_transform(broomfield_precincts, crs=prg)) %>%
rbind(st_transform(westminster_city, crs=prg))
# Marshall fire boundary
wfigs_fire = st_read("../GIS_inputs_destruction_fireboundary/WFIGS_-_Wildland_Fire_Perimeters_Full_History/FH_Perimeter.shp")
## Reading layer `FH_Perimeter' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\WFIGS_-_Wildland_Fire_Perimeters_Full_History\FH_Perimeter.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 9500 features and 108 fields (with 4 geometries empty)
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -165.8152 ymin: 13.38021 xmax: 144.6906 ymax: 69.05174
## Geodetic CRS: WGS 84
# filter fire to marshall fire only; update crs
(marshall_fire = wfigs_fire %>% filter(poly_Incid == "Marshall") %>% st_transform(crs = st_crs(prg)))
## Simple feature collection with 1 feature and 108 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -105.2333 ymin: 39.92924 xmax: -105.1314 ymax: 39.98638
## CRS: +proj=longlat +datum=WGS84 +no_defs
## OBJECTID poly_Incid poly_Featu poly_MapMe poly_GISAc
## 1 14624 Marshall Wildfire Daily Fire Perimeter Mixed Methods 6079.904
## poly_Creat poly_DateC poly_Polyg poly_Acres
## 1 2022-01-03 2022-01-06 <NA> 6079.882
## poly_Globa poly_Sourc irwin_ABCD irwin_ADSP
## 1 {50133057-91F9-4FC7-92AA-34A69C58684C} 2022 NIFS <NA> FIREREPORTING
## irwin_Calc irwin_Cont irwin_Co_1 irwin_Dail irwin_Disc irwin_Disp irwin_Esti
## 1 NA 2022-01-03 2022-01-05 6026 1 COFTC 2349076
## irwin_FFRe irwin_FF_1 irwin_FF_2 irwin_Fire irwin_Fi_1 irwin_Fi_2 irwin_Fi_3
## 1 <NA> <NA> <NA> Minimal Smoldering <NA> <NA>
## irwin_Fi_4 irwin_Fi_5 irwin_Fi_6 irwin_Fi_7 irwin_Fi_8 irwin_Fi_9
## 1 Undetermined <NA> <NA> PER9 <NA> 2021-12-30
## irwin_Fi10 irwin_Fi11 irwin_FSCo irwin_FSFu irwin_FSMo irwin_FSPo
## 1 Type 1 Incident 2022-01-08 0 100 0 0
## irwin_FSJo irwin_FSOv irwin_GACC irwin_ICS2 irwin_IC_1 irwin_IC_2 irwin_IC_3
## 1 PN 1522 RMCC 2022-01-07 2022-01-07 2022-01-07 F
## irwin_Inci irwin_In_1 irwin_In_2 irwin_In_3
## 1 Type 1 Team Marshall 5 Miles Southeast of Boulder Colorado WF
## irwin_In_4 irwin_Init irwin_In_5 irwin_In_6 irwin_In_7
## 1 FI 39.93435 -105.1362 6000 <NA>
## irwin_Irwi irwin_IsFi irwin_Is_1 irwin_IsFS
## 1 {C63FC371-BC70-4615-841B-B0838C21064F} NA 0 1
## irwin_IsMu irwin_IsRe irwin_IsTr irwin_IsUn irwin_Loca irwin_Perc irwin_Pe_1
## 1 0 0 0 0 000995 100 100
## irwin_POOC irwin_PO_1 irwin_POOD irwin_POOF irwin_POOJ irwin_PO_2 irwin_PO_3
## 1 <NA> Boulder COFTC 08013 <NA> <NA> <NA>
## irwin_POOL irwin_PO_4 irwin_PO_5 irwin_PO_6 irwin_PO_7 irwin_PO_8 irwin_PO_9
## 1 Private Private <NA> <NA> <NA> <NA> NA
## irwin_PO10 irwin_POOP irwin_PO11 irwin_PO12 irwin_POOS irwin_Pred irwin_Pr_1
## 1 <NA> RM12 C&L COBLX US-CO <NA> <NA>
## irwin_Prim irwin_Seco irwin_Tota
## 1 Short Grass (1 foot) Timber (Grass and Understory) 57
## irwin_Uniq irwin_WFDS irwin_Crea irwin_Modi irwin_IsDi irwin_Orga
## 1 2021-COBLX-000995 No Decision wildcad wildcad 0 <NA>
## irwin_Stra irwin_Glob irwin_Sour irwin_Arch irwin_Mo_1 irwin_Cr_1
## 1 <NA> <NA> IRWIN NA 2022-02-17 2021-12-30
## GlobalID irwin_IsCp irwin_CpxN irwin_CpxI
## 1 {D2E42D65-AA00-4926-9B3A-09C8620C7A78} 0 <NA> <NA>
## SHAPE_Leng SHAPE_Area geometry
## 1 0.5953818 -0.002593346 MULTIPOLYGON (((-105.1558 3...
# BUILDINGS
# Read in destroyed home and building sites
destroyed_homes = st_read("../GIS_inputs_destruction_fireboundary/output_damage_files/destroyed_homes.shp")
## Reading layer `destroyed_homes' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\output_damage_files\destroyed_homes.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 1030 features and 25 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.23 ymin: 39.93075 xmax: -105.135 ymax: 40.01238
## Geodetic CRS: WGS 84
damaged_homes = st_read("../GIS_inputs_destruction_fireboundary/output_damage_files/damaged_homes.shp")
## Reading layer `damaged_homes' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\output_damage_files\damaged_homes.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 143 features and 25 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.2309 ymin: 39.93081 xmax: -105.1439 ymax: 40.01083
## Geodetic CRS: WGS 84
destroyed_businesses = st_read("../GIS_inputs_destruction_fireboundary/output_damage_files/destroyed_businesses.shp")
## Reading layer `destroyed_businesses' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\output_damage_files\destroyed_businesses.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 7 features and 24 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.1651 ymin: 39.94129 xmax: -105.1384 ymax: 39.95788
## Geodetic CRS: WGS 84
damaged_businesses = st_read("../GIS_inputs_destruction_fireboundary/output_damage_files/damaged_businesses.shp")
## Reading layer `damaged_businesses' from data source
## `C:\Users\zclem\Documents\GitHub\marshall_fires\GIS_inputs_destruction_fireboundary\output_damage_files\damaged_businesses.shp'
## using driver `ESRI Shapefile'
## Simple feature collection with 27 features and 24 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.2132 ymin: 39.94129 xmax: -105.1384 ymax: 40.01153
## Geodetic CRS: WGS 84
# create 'type' column; select only the columns needed; make projection string match
(destroyed_homes = destroyed_homes %>% mutate(type = "destroyed - residential") %>% dplyr::rename(jurisdiction = JURISDICTI) %>% dplyr::select(jurisdiction, type, latlong, geometry) %>% st_transform(crs = st_crs(prg)))
## Simple feature collection with 1030 features and 3 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.23 ymin: 39.93075 xmax: -105.135 ymax: 40.01238
## CRS: +proj=longlat +datum=WGS84 +no_defs
## First 10 features:
## jurisdiction type latlong
## 1 Louisville destroyed - residential 39.977609,-105.163733
## 2 Louisville destroyed - residential 39.97682,-105.163673
## 3 Louisville destroyed - residential 39.97682,-105.163673
## 4 Louisville destroyed - residential 39.97682,-105.163673
## 5 Louisville destroyed - residential 39.9777,-105.163555
## 6 Louisville destroyed - residential 39.97771569708369,-105.16356686407488
## 7 Louisville destroyed - residential 39.9777426,-105.16414479323862
## 8 Louisville destroyed - residential 39.97769007192905,-105.16382127167078
## 9 Louisville destroyed - residential 39.97774709125106,-105.16359059222464
## 10 Louisville destroyed - residential 39.976826,-105.162003
## geometry
## 1 POINT (-105.1637 39.97761)
## 2 POINT (-105.1637 39.97682)
## 3 POINT (-105.1637 39.97682)
## 4 POINT (-105.1637 39.97682)
## 5 POINT (-105.1636 39.9777)
## 6 POINT (-105.1636 39.97772)
## 7 POINT (-105.1641 39.97774)
## 8 POINT (-105.1638 39.97769)
## 9 POINT (-105.1636 39.97775)
## 10 POINT (-105.162 39.97683)
(destroyed_businesses = destroyed_businesses %>% mutate(type = "destroyed - non-residential") %>% dplyr::rename(jurisdiction = Jurisdicti) %>% dplyr::select(jurisdiction, type, latlong, geometry) %>% st_transform(crs = st_crs(prg)))
## Simple feature collection with 7 features and 3 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.1651 ymin: 39.94129 xmax: -105.1384 ymax: 39.95788
## CRS: +proj=longlat +datum=WGS84 +no_defs
## jurisdiction type latlong
## 1 Louisville destroyed - non-residential 39.9578801,-105.1384268
## 2 Louisville destroyed - non-residential 39.9578801,-105.1384268
## 3 Louisville destroyed - non-residential 39.9578801,-105.1384268
## 4 Louisville destroyed - non-residential 39.9578801,-105.1384268
## 5 Superior destroyed - non-residential 39.9412918,-105.1511702
## 6 Superior destroyed - non-residential 39.95469525,-105.16507092452144
## 7 Superior destroyed - non-residential 39.9552909,-105.1639793
## geometry
## 1 POINT (-105.1384 39.95788)
## 2 POINT (-105.1384 39.95788)
## 3 POINT (-105.1384 39.95788)
## 4 POINT (-105.1384 39.95788)
## 5 POINT (-105.1512 39.94129)
## 6 POINT (-105.1651 39.9547)
## 7 POINT (-105.164 39.95529)
(damaged_homes = damaged_homes %>% mutate(type = "damaged - residential") %>% dplyr::rename(jurisdiction = JURISDICTI) %>% dplyr::select(jurisdiction, type, latlong, geometry) %>% st_transform(crs = st_crs(prg)))
## Simple feature collection with 143 features and 3 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.2309 ymin: 39.93081 xmax: -105.1439 ymax: 40.01083
## CRS: +proj=longlat +datum=WGS84 +no_defs
## First 10 features:
## jurisdiction type latlong
## 1 Louisville damaged - residential 39.97773139416737,-105.16357872814976
## 2 Louisville damaged - residential 39.97676665,-105.16478214953415
## 3 Louisville damaged - residential 39.97146875,-105.15768883488613
## 4 Louisville damaged - residential 39.9709326700897,-105.15778031176009
## 5 Louisville damaged - residential 39.970950945310186,-105.15756361107347
## 6 Louisville damaged - residential 39.97375085,-105.15190285613079
## 7 Louisville damaged - residential 39.97774621133401,-105.16066240508351
## 8 Louisville damaged - residential 39.9713432,-105.1663380652309
## 9 Louisville damaged - residential 39.978311903473696,-105.15760904854133
## 10 Louisville damaged - residential 39.97795856454827,-105.15977564051832
## geometry
## 1 POINT (-105.1636 39.97773)
## 2 POINT (-105.1648 39.97677)
## 3 POINT (-105.1577 39.97147)
## 4 POINT (-105.1578 39.97093)
## 5 POINT (-105.1576 39.97095)
## 6 POINT (-105.1519 39.97375)
## 7 POINT (-105.1607 39.97775)
## 8 POINT (-105.1663 39.97134)
## 9 POINT (-105.1576 39.97831)
## 10 POINT (-105.1598 39.97796)
(damaged_business = damaged_businesses %>% mutate(type = "damaged - non-residential") %>% dplyr::rename(jurisdiction = Jurisdicti) %>% dplyr::select(jurisdiction, type, latlong, geometry) %>% st_transform(crs = st_crs(prg)))
## Simple feature collection with 27 features and 3 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.2132 ymin: 39.94129 xmax: -105.1384 ymax: 40.01153
## CRS: +proj=longlat +datum=WGS84 +no_defs
## First 10 features:
## jurisdiction type
## 1 Louisville damaged - non-residential
## 2 Louisville damaged - non-residential
## 3 Louisville damaged - non-residential
## 4 Louisville damaged - non-residential
## 5 Louisville damaged - non-residential
## 6 Louisville damaged - non-residential
## 7 Louisville damaged - non-residential
## 8 Louisville damaged - non-residential
## 9 Louisville damaged - non-residential
## 10 Louisville damaged - non-residential
## latlong geometry
## 1 39.9649982,-105.16525087536996 POINT (-105.1653 39.965)
## 2 39.96509899678052,-105.16649768361437 POINT (-105.1665 39.9651)
## 3 39.9706,-105.166592 POINT (-105.1666 39.9706)
## 4 39.96818493380565,-105.16767989691307 POINT (-105.1677 39.96818)
## 5 39.975743702971805,-105.16600406830685 POINT (-105.166 39.97574)
## 6 39.964466787878784,-105.16362403030303 POINT (-105.1636 39.96447)
## 7 39.9578801,-105.1384268 POINT (-105.1384 39.95788)
## 8 39.9578801,-105.1384268 POINT (-105.1384 39.95788)
## 9 39.9578801,-105.1384268 POINT (-105.1384 39.95788)
## 10 39.9578801,-105.1384268 POINT (-105.1384 39.95788)
# combine into one df, using the "," indicator keeps decimal points
(destroyed_damaged =rbind(destroyed_homes,destroyed_businesses,damaged_homes,damaged_business) %>% separate(latlong, c("lat","long"), ",", convert = FALSE))
## Simple feature collection with 1207 features and 4 fields
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -105.2309 ymin: 39.93075 xmax: -105.135 ymax: 40.01238
## CRS: +proj=longlat +datum=WGS84 +no_defs
## First 10 features:
## jurisdiction type lat long
## 1 Louisville destroyed - residential 39.977609 -105.163733
## 2 Louisville destroyed - residential 39.97682 -105.163673
## 3 Louisville destroyed - residential 39.97682 -105.163673
## 4 Louisville destroyed - residential 39.97682 -105.163673
## 5 Louisville destroyed - residential 39.9777 -105.163555
## 6 Louisville destroyed - residential 39.97771569708369 -105.16356686407488
## 7 Louisville destroyed - residential 39.9777426 -105.16414479323862
## 8 Louisville destroyed - residential 39.97769007192905 -105.16382127167078
## 9 Louisville destroyed - residential 39.97774709125106 -105.16359059222464
## 10 Louisville destroyed - residential 39.976826 -105.162003
## geometry
## 1 POINT (-105.1637 39.97761)
## 2 POINT (-105.1637 39.97682)
## 3 POINT (-105.1637 39.97682)
## 4 POINT (-105.1637 39.97682)
## 5 POINT (-105.1636 39.9777)
## 6 POINT (-105.1636 39.97772)
## 7 POINT (-105.1641 39.97774)
## 8 POINT (-105.1638 39.97769)
## 9 POINT (-105.1636 39.97775)
## 10 POINT (-105.162 39.97683)
# change lat and long to floats
destroyed_damaged$lat = as.double(destroyed_damaged$lat)
destroyed_damaged$long = as.double(destroyed_damaged$long)
# create "location status" column for mapping]
sensors$loc_status = str_c(sensors$Location, " - ", sensors$Status)
# turned this on to get code to run -- should update classification code in the extraction file
sensors = na.omit(sensors)
# set factors to false
options(stringsAsFactors = FALSE)
# create pallete for destroyed and damaged building colors for leaflet map
pal = colorFactor(c("red","orange","magenta","light pink"), domain = unique(destroyed_damaged$type))
## create indoor/outdoor icons
sensorIcon = awesomeIconList(
"inside" = makeAwesomeIcon(
icon = "home",
iconColor = "white",
markerColor = "blue",
library = "fa"
),
"outside" = makeAwesomeIcon(
icon = "tree",
iconColor = "white",
markerColor = "green",
library = "fa"
))
sensors
## X ID fire_period post_fire_period pre_fire_period
## 1 1 3062 100.000000 100.00000 100.00000
## 2 2 4583 100.000000 99.93082 100.00000
## 3 3 4643 100.000000 99.33126 100.00000
## 4 4 15663 100.000000 99.52727 100.00000
## 5 5 27120 48.333333 99.86740 100.00000
## 6 6 30743 74.000000 73.96770 100.00000
## 7 7 33109 53.333333 99.70022 100.00000
## 8 8 38623 7.333333 95.18621 100.00000
## 10 10 45325 99.666667 88.22783 100.00000
## 11 11 45369 99.666667 99.94235 100.00000
## 13 13 50133 100.000000 100.00000 100.00000
## 14 14 53729 100.000000 99.84434 100.00000
## 16 16 56209 100.000000 94.93255 96.66667
## 17 17 56431 86.666667 99.87893 100.00000
## 18 18 60507 100.000000 99.93188 100.00000
## 19 19 60517 100.000000 99.82705 100.00000
## 20 20 64671 100.000000 99.81552 100.00000
## 21 21 64967 22.333333 99.86740 100.00000
## 22 22 65669 17.666667 99.89046 100.00000
## 23 23 66173 100.000000 99.96541 100.00000
## 24 24 67917 99.333333 99.96541 100.00000
## 25 25 68375 7.333333 99.25055 100.00000
## 26 26 69103 4.000000 99.72328 100.00000
## 28 28 70241 100.000000 99.96541 100.00000
## 29 29 70911 100.000000 94.94408 100.00000
## 31 31 72931 100.000000 99.94811 100.00000
## 34 34 79121 7.333333 92.70725 100.00000
## 35 35 79907 38.000000 99.99423 100.00000
## 36 36 81149 100.000000 99.59068 100.00000
## 37 37 82879 100.000000 99.68292 100.00000
## 38 38 83423 100.000000 98.81817 100.00000
## 39 39 83597 97.333333 99.58492 100.00000
## 40 40 83661 100.000000 95.07091 100.00000
## 41 41 83785 4.000000 99.72328 100.00000
## 42 42 84169 100.000000 99.91929 100.00000
## 43 43 84585 76.000000 98.76780 100.00000
## 44 44 84697 1.333333 96.62170 100.00000
## 45 45 84777 10.666667 63.97674 100.00000
## 46 46 85643 100.000000 100.00000 100.00000
## 47 47 85999 14.000000 97.89000 100.00000
## 48 48 86477 13.666667 98.16672 100.00000
## 49 49 86637 99.666667 99.78669 100.00000
## 50 50 87465 100.000000 99.98847 100.00000
## 51 51 87531 100.000000 99.38314 100.00000
## 52 52 88035 100.000000 99.52150 100.00000
## 53 53 88533 100.000000 99.97694 100.00000
## 55 55 90255 3.666667 99.64257 100.00000
## 56 56 90297 89.000000 90.30900 100.00000
## 57 57 91877 100.000000 99.94811 100.00000
## 58 58 92321 100.000000 99.97117 100.00000
## 59 59 94411 100.000000 88.86199 100.00000
## 60 60 95591 100.000000 100.00000 100.00000
## 61 61 95831 99.000000 87.35155 100.00000
## 62 62 96731 100.000000 100.00000 100.00000
## 63 63 97647 100.000000 99.93658 100.00000
## 64 64 100521 100.000000 99.60798 100.00000
## 65 65 101003 4.000000 99.40044 100.00000
## 66 66 102614 100.000000 99.47538 100.00000
## 67 67 102724 100.000000 99.82128 100.00000
## 69 69 103724 7.333333 98.20708 100.00000
## 70 70 103820 49.000000 46.41992 100.00000
## 71 71 104594 99.666667 99.78093 100.00000
## 72 72 106214 22.333333 99.85587 100.00000
## 73 73 106492 48.666667 98.85276 100.00000
## 74 74 106800 100.000000 99.67139 100.00000
## 75 75 106906 99.666667 99.99423 100.00000
## 76 76 107466 63.000000 99.96541 100.00000
## 78 78 108558 99.666667 99.77516 100.00000
## 80 80 110042 100.000000 100.00000 100.00000
## 81 81 110122 20.666667 99.96541 100.00000
## 82 82 111834 56.666667 99.34855 100.00000
## 83 83 112082 56.666667 99.59068 100.00000
## 84 84 112522 100.000000 97.51528 100.00000
## 85 85 112606 99.666667 97.58446 100.00000
## 86 86 112974 100.000000 99.98847 100.00000
## 87 87 113054 100.000000 99.96541 100.00000
## 88 88 113266 100.000000 99.87317 100.00000
## 93 93 115303 99.666667 99.94811 100.00000
## 94 94 117647 100.000000 99.82705 100.00000
## 95 95 119547 99.666667 99.44656 100.00000
## 97 97 120477 100.000000 99.96541 100.00000
## 99 99 120847 53.333333 99.98812 100.00000
## 101 101 121197 100.000000 99.98656 100.00000
## 103 103 122071 98.333333 99.00842 100.00000
## 106 106 122585 4.000000 96.54675 100.00000
## 115 115 124319 53.333333 99.46962 100.00000
## 200 200 141956 45.000000 64.57397 100.00000
## Status Month
## 1 Complete data throughout fire period Before or during 12-2021
## 2 Complete data throughout fire period Before or during 12-2021
## 3 Complete data throughout fire period Before or during 12-2021
## 4 Complete data throughout fire period Before or during 12-2021
## 5 Sensor offline during fire, returned online Before or during 12-2021
## 6 Sensor offline during fire, did not return online Before or during 12-2021
## 7 Sensor offline during fire, returned online Before or during 12-2021
## 8 Sensor offline during fire, returned online Before or during 12-2021
## 10 Complete data throughout fire period Before or during 12-2021
## 11 Complete data throughout fire period Before or during 12-2021
## 13 Complete data throughout fire period Before or during 12-2021
## 14 Complete data throughout fire period Before or during 12-2021
## 16 Complete data throughout fire period Before or during 12-2021
## 17 Complete data throughout fire period Before or during 12-2021
## 18 Complete data throughout fire period Before or during 12-2021
## 19 Complete data throughout fire period Before or during 12-2021
## 20 Complete data throughout fire period Before or during 12-2021
## 21 Sensor offline during fire, returned online Before or during 12-2021
## 22 Sensor offline during fire, returned online Before or during 12-2021
## 23 Complete data throughout fire period Before or during 12-2021
## 24 Complete data throughout fire period Before or during 12-2021
## 25 Sensor offline during fire, returned online Before or during 12-2021
## 26 Sensor offline during fire, returned online Before or during 12-2021
## 28 Complete data throughout fire period Before or during 12-2021
## 29 Complete data throughout fire period Before or during 12-2021
## 31 Complete data throughout fire period Before or during 12-2021
## 34 Sensor offline during fire, returned online Before or during 12-2021
## 35 Sensor offline during fire, returned online Before or during 12-2021
## 36 Complete data throughout fire period Before or during 12-2021
## 37 Complete data throughout fire period Before or during 12-2021
## 38 Complete data throughout fire period Before or during 12-2021
## 39 Complete data throughout fire period Before or during 12-2021
## 40 Complete data throughout fire period Before or during 12-2021
## 41 Sensor offline during fire, returned online Before or during 12-2021
## 42 Complete data throughout fire period Before or during 12-2021
## 43 Complete data throughout fire period Before or during 12-2021
## 44 Sensor offline during fire, returned online Before or during 12-2021
## 45 Sensor offline during fire, did not return online Before or during 12-2021
## 46 Complete data throughout fire period Before or during 12-2021
## 47 Sensor offline during fire, returned online Before or during 12-2021
## 48 Sensor offline during fire, returned online Before or during 12-2021
## 49 Complete data throughout fire period Before or during 12-2021
## 50 Complete data throughout fire period Before or during 12-2021
## 51 Complete data throughout fire period Before or during 12-2021
## 52 Complete data throughout fire period Before or during 12-2021
## 53 Complete data throughout fire period Before or during 12-2021
## 55 Sensor offline during fire, returned online Before or during 12-2021
## 56 Complete data throughout fire period Before or during 12-2021
## 57 Complete data throughout fire period Before or during 12-2021
## 58 Complete data throughout fire period Before or during 12-2021
## 59 Complete data throughout fire period Before or during 12-2021
## 60 Complete data throughout fire period Before or during 12-2021
## 61 Complete data throughout fire period Before or during 12-2021
## 62 Complete data throughout fire period Before or during 12-2021
## 63 Complete data throughout fire period Before or during 12-2021
## 64 Complete data throughout fire period Before or during 12-2021
## 65 Sensor offline during fire, returned online Before or during 12-2021
## 66 Complete data throughout fire period Before or during 12-2021
## 67 Complete data throughout fire period Before or during 12-2021
## 69 Sensor offline during fire, returned online Before or during 12-2021
## 70 Sensor offline during fire, did not return online Before or during 12-2021
## 71 Complete data throughout fire period Before or during 12-2021
## 72 Sensor offline during fire, returned online Before or during 12-2021
## 73 Sensor offline during fire, returned online Before or during 12-2021
## 74 Complete data throughout fire period Before or during 12-2021
## 75 Complete data throughout fire period Before or during 12-2021
## 76 Sensor offline during fire, returned online Before or during 12-2021
## 78 Complete data throughout fire period Before or during 12-2021
## 80 Complete data throughout fire period Before or during 12-2021
## 81 Sensor offline during fire, returned online Before or during 12-2021
## 82 Sensor offline during fire, returned online Before or during 12-2021
## 83 Sensor offline during fire, returned online Before or during 12-2021
## 84 Complete data throughout fire period Before or during 12-2021
## 85 Complete data throughout fire period Before or during 12-2021
## 86 Complete data throughout fire period Before or during 12-2021
## 87 Complete data throughout fire period Before or during 12-2021
## 88 Complete data throughout fire period Before or during 12-2021
## 93 Complete data throughout fire period Before or during 12-2021
## 94 Complete data throughout fire period Before or during 12-2021
## 95 Complete data throughout fire period Before or during 12-2021
## 97 Complete data throughout fire period Before or during 12-2021
## 99 Sensor offline during fire, returned online Before or during 12-2021
## 101 Complete data throughout fire period Before or during 12-2021
## 103 Complete data throughout fire period Before or during 12-2021
## 106 Sensor offline during fire, returned online Before or during 12-2021
## 115 Sensor offline during fire, returned online Before or during 12-2021
## 200 Sensor offline during fire, did not return online Before or during 12-2021
## Lon Lat Name Location
## 1 -105.2713 40.04077 Kalmia outside
## 2 -105.1680 40.13818 Clover Basin Consulting outside
## 3 -105.2515 39.98943 JENSA outside
## 4 -105.1550 40.10230 FRCC Niwot location outside
## 5 -105.2600 39.98123 Table Mesa Court outside
## 6 -105.1283 40.13157 Creekside outside
## 7 -105.2622 39.98305 Bear outside
## 8 -105.1749 39.98661 76th & S Bldr Rd outside
## 10 -105.1647 40.06149 Kincross Way outside
## 11 -105.1649 40.06169 Kincross Way inside
## 13 -105.1727 40.10403 HRNS-TEST inside
## 14 -105.0457 39.99003 Anthem Ranch outside
## 16 -105.2691 39.97224 Stony Hill Road outside
## 17 -105.2821 40.03435 Grape_and_Broadway outside
## 18 -105.2694 40.01557 Bella outside
## 19 -105.1199 39.87919 STANDLEY LAKE outside
## 20 -105.2709 40.04019 K2055i inside
## 21 -105.2458 39.98447 Majestic Heights, Boulder, CO outside
## 22 -105.1477 39.97686 Louisville outside
## 23 -105.1214 39.94361 Terracina Broomfield Apartments outside
## 24 -105.2833 40.00194 Chautauqua Heights outside
## 25 -105.1500 39.98119 Louisville Air Sensor outside
## 26 -105.3219 40.04352 Hawk Lane outside
## 28 -105.2804 40.02177 MKHOME2 outside
## 29 -105.2508 39.97368 Shanahan Ridge inside
## 31 -105.2279 39.99787 AZTEC CENTRAL outside
## 34 -105.1741 39.98653 Kitchen inside
## 35 -105.2716 40.03448 Vista Drive inside
## 36 -105.1549 39.97095 ntsky outside
## 37 -105.2401 39.98775 Mango Manor Coop inside
## 38 -105.2849 40.05431 Shining Mountain Waldorf School outside
## 39 -105.2563 39.97676 South Boulder outside
## 40 -105.1233 40.13444 Creekside and Sunset outside
## 41 -105.3218 40.04347 Hawk Lane inside
## 42 -105.2312 40.01628 Snowcave inside
## 43 -105.2308 39.99910 645 Manhattan Place inside
## 44 -105.2798 39.93268 Eldorado Springs outside
## 45 -105.1660 39.95349 Superior Town Hall outside
## 46 -105.3186 40.02700 Lower Seven Hills outside
## 47 -105.1554 39.93897 Superior - North Pool outside
## 48 -105.1552 39.92155 Superior - South Pool outside
## 49 -105.2515 39.97368 Viele Lake outside
## 50 -105.2845 40.06458 Dakota Ridge outside
## 51 -105.2394 39.98681 The Haunted Pickle Mansion inside
## 52 -105.1928 40.07116 Gunbarrel Green Spotted Horse (inside) inside
## 53 -105.1489 39.99006 Purple house indoors inside
## 55 -105.2551 39.98033 Miami Way outside
## 56 -105.0926 39.96745 South Pointe outside
## 57 -105.1451 39.98836 Eisenhower Dr outside
## 58 -105.1491 39.99013 Purple House outside
## 59 -105.2689 40.03923 Parkside outside
## 60 -105.2893 40.05224 Wonderland Lake outside
## 61 -105.1641 40.06173 Gunbarrel Hill outside
## 62 -105.2875 40.02247 Concord Ave outside
## 63 -105.0223 39.92575 Broomfield - Calkins Pl outside
## 64 -105.1932 40.07457 PurpleHome inside
## 65 -105.3151 40.04732 House inside
## 66 -105.0764 40.04520 Ponderosa outside
## 67 -105.1251 40.02127 Blue Heron Estates outside
## 69 -105.1499 39.97538 Home inside
## 70 -105.2813 40.00443 Casper’s House inside
## 71 -105.2007 40.07604 Bean Mountain Ln outside
## 72 -105.2458 39.98451 Majestic Heights inside
## 73 -105.3364 40.07924 Boulder Heights inside
## 74 -105.1249 40.02118 Blue Heron Estates inside
## 75 -105.2670 39.97523 1333 Wildwood Ct., Boulder, CO 80305 Inside inside
## 76 -105.2559 39.98311 Yale Road Home inside
## 78 -105.1720 40.09753 CountryCreek2 outside
## 80 -105.2562 40.03099 Valmont and 29th outside
## 81 -105.2427 39.98229 Hanover Ave. inside
## 82 -105.2451 40.09622 Snead Ct outside
## 83 -105.2549 40.09183 Eagle Ct. outside
## 84 -105.0474 39.82646 Linda Park Addition outside
## 85 -105.2670 39.97515 1333 Wildwood Ct., Boulder, CO 80305 outside
## 86 -105.1615 40.13421 Summerlin_Lane outside
## 87 -105.2837 40.04738 Quince Ave inside
## 88 -105.2803 40.04895 Redwood & Riverside outside
## 93 -105.0864 39.90192 Amli Arista outside
## 94 -105.1197 39.87945 Stanley Lake inside
## 95 -105.0439 39.95181 Broadlands outside
## 97 -105.1246 39.99547 Bohannan house outside
## 99 -105.2659 39.97925 Yard outside
## 101 -105.2376 40.00849 Hancock Drive outside
## 103 -105.1962 40.04964 EMSonJay outside
## 106 -105.2556 39.98821 Bates inside
## 115 -105.2567 39.98218 Hartford Drive inside
## 200 -105.1550 39.97581 Nighthawk Resident inside
## address
## 1 2065, Kalmia Avenue, Boulder, Boulder County, Colorado, 80304, United States
## 2 5017, William Place, Longmont, Boulder County, Colorado, 80503, United States
## 3 3565, Eastman Avenue, Boulder, Boulder County, Colorado, 80305, United States
## 4 7037, Longview Drive, Boulder, Niwot, Boulder County, Colorado, 80503, United States
## 5 2799, Table Mesa Court, Boulder, Boulder County, Colorado, 80305, United States
## 6 2443, Eagleview Circle, Longmont, Boulder County, Colorado, 80504, United States
## 7 732, Ithaca Drive, Boulder, Boulder County, Colorado, 80305, United States
## 8 7618, South Boulder Road, Boulder, Boulder County, Colorado, 80303, United States
## 10 8059, Kincross Way, Heatherwood, Boulder, Boulder County, Colorado, 80301, United States
## 11 8059, Kincross Way, Heatherwood, Boulder, Boulder County, Colorado, 80301, United States
## 13 My Mama's Pie, Murray Street, Boulder, Niwot, Boulder County, Colorado, 80544, United States
## 14 4512, Silver Mountain Loop, Anthem, Broomfield, Colorado, 80032, United States
## 16 1916, Stony Hill Road, Boulder, Boulder County, Colorado, 80305, United States
## 17 Broadway & Grape Ave, Broadway, Washington Village, Boulder, Boulder County, Colorado, 80306, United States
## 18 1933, Grove Street, Boulder, Boulder County, Colorado, 80302, United States
## 19 10098, Oak Court, Walnut Grove, Westminster, Jefferson County, Colorado, 80021, United States
## 20 2099, Kalmia Avenue, Boulder, Boulder County, Colorado, 80304, United States
## 21 605, South 42nd Street, Boulder, Boulder County, Colorado, 80305, United States
## 22 661, Cleveland Avenue, Louisville, Boulder County, Colorado, 80027, United States
## 23 Building 30, 13630, Via Varra, Overlook District, Broomfield, Colorado, 80020, United States
## 24 847, Cascade Avenue, Central Boulder - University Hill, Boulder, Boulder County, Colorado, 80302, United States
## 25 599, West Sagebrush Court, Louisville, Boulder County, Colorado, 80027, United States
## 26 67, Hawk Lane, Boulder County, Colorado, 80304, United States
## 28 2357, 13th Street, Washington Village, Boulder, Boulder County, Colorado, 80304, United States
## 29 Boulder, Boulder County, Colorado, 80305, United States
## 31 556, Aztec Drive, Boulder, Boulder County, Colorado, 80303, United States
## 34 South Boulder Road, Boulder, Boulder County, Colorado, 80027, United States
## 35 1970, Vista Drive, Boulder, Boulder County, Colorado, 80304, United States
## 36 775, West Lois Court, Louisville, Boulder County, Colorado, 80027, United States
## 37 4662, Ingram Court, Boulder, Boulder County, Colorado, 80305, United States
## 38 Shining Mountain Waldorf School, 999, Violet Avenue, Boulder, Boulder County, Colorado, 80304, United States
## 39 1264, Ithaca Drive, Boulder, Boulder County, Colorado, 80305, United States
## 40 1839, Caleta Trail, Longmont, Boulder County, Colorado, 80504, United States
## 41 67, Hawk Lane, Boulder County, Colorado, 80304, United States
## 42 1707, Commerce Street, Boulder, Boulder County, Colorado, 80301, United States
## 43 645, Manhattan Drive, Boulder, Boulder County, Colorado, 80303, United States
## 44 Eldorado Springs Pool, Artesian Drive, Eldorado Springs, Boulder County, Colorado, 80025, United States
## 45 Superior Municipal Court, East William Street, Original Superior, Superior, Boulder County, Colorado, 80027, United States
## 46 442, Seven Hills Drive, Seven Hills, Boulder County, Colorado, 80304, United States
## 47 Town of Superior North Pool Facility, South Indiana Street, Rock Creek Ranch, Boulder, Superior, Boulder County, Colorado, 80027, United States
## 48 Superior South Pool, Huron Peak Drive, Rock Creek Ranch II, Boulder, Boulder County, Colorado, 80027, United States
## 49 3270, Everett Drive, Boulder, Boulder County, Colorado, 80305, United States
## 50 4949, Broadway, Boulder County, Colorado, 80304, United States
## 51 4682, Ingram Court, Boulder, Boulder County, Colorado, 80305, United States
## 52 5310, Spotted Horse Trail, Boulder, Boulder County, Colorado, 80301, United States
## 53 1998, Eisenhower Drive, Louisville, Boulder County, Colorado, 80027, United States
## 55 950, Miami Way, Boulder, Boulder County, Colorado, 80305, United States
## 56 South Point Drive, South Pointe, Lafayette, Boulder County, Colorado, 80026-2872, United States
## 57 370, Eisenhower Drive, Louisville, Boulder County, Colorado, 80027, United States
## 58 1998, Eisenhower Drive, Louisville, Boulder County, Colorado, 80027, United States
## 59 22nd Street, Boulder, Boulder County, Colorado, 80304, United States
## 60 501, Utica Avenue, Boulder, Boulder County, Colorado, 80303, United States
## 61 8067, Kincross Way, Heatherwood, Boulder, Boulder County, Colorado, 80301, United States
## 62 Concord Alley, Washington Village, Boulder, Boulder County, Colorado, 80302, United States
## 63 2830, Calkins Place, Broomfield, Colorado, 80020, United States
## 64 5583, Homestead Way, Boulder, Boulder County, Colorado, 80301, United States
## 65 230, Timber Lane, Pine Brook Hill, Boulder County, Colorado, 80304, United States
## 66 11821, Flatiron Drive, Leyner, Erie, Boulder County, Colorado, 80026, United States
## 67 2432, Ginny Way, Blue Heron Estates, Lafayette, Boulder County, Colorado, 80026, United States
## 69 681, Spruce Circle, Louisville, Boulder County, Colorado, 80027, United States
## 70 941, Lincoln Place, Central Boulder - University Hill, Boulder, Boulder County, Colorado, 80302, United States
## 71 6655, Bean Mountain Lane, Boulder, Boulder County, Colorado, 80301, United States
## 72 605, South 42nd Street, Boulder, Boulder County, Colorado, 80305, United States
## 73 253, Deer Trail Road, Lazy Acres, Boulder County, Colorado, 80455, United States
## 74 2430, Ginny Way, Blue Heron Estates, Lafayette, Boulder County, Colorado, 80026, United States
## 75 1313, Wildwood Court, Boulder, Boulder County, Colorado, 80305, United States
## 76 720, Yale Road, Boulder, Boulder County, Colorado, 80305, United States
## 78 7777, Country Creek Drive, Boulder, Niwot, Boulder County, Colorado, 80503, United States
## 80 3185, 29th Street, Boulder, Boulder County, Colorado, 80301, United States
## 81 770, South 45th Street, Boulder, Boulder County, Colorado, 80305, United States
## 82 6732, Snead Court, Boulder County, Colorado, 80503, United States
## 83 6430, Eagle Court, Boulder County, Colorado, 80503, United States
## 84 7157, Winona Court, Westminster, Adams County, Colorado, 80030, United States
## 85 1313, Wildwood Court, Boulder, Boulder County, Colorado, 80305, United States
## 86 Longmont, Boulder County, Colorado, 80503, United States
## 87 1052, Quince Avenue, Boulder, Boulder County, Colorado, 80304, United States
## 88 1355, Redwood Avenue, Boulder, Boulder County, Colorado, 80304, United States
## 93 7971, Uptown Avenue, Arista, Broomfield, Colorado, 80021, United States
## 94 10081, Oak Street, Westminster, Jefferson County, Colorado, 80021, United States
## 95 4342, Augusta Drive, Broadlands, Broomfield, Colorado, 80023, United States
## 97 384 Driftwood Circle, 384, Driftwood Circle, Waneka Landing, Lafayette, Boulder County, Colorado, 80026, United States
## 99 Magic Mansion, 2267, Holyoke Drive, Boulder, Boulder County, Colorado, 80305, United States
## 101 3878, Hancock Drive, Boulder, Boulder County, Colorado, 80303, United States
## 103 Jay Road, Boulder, Boulder County, Colorado, 80301, United States
## 106 350, Bates Avenue, Boulder, Boulder County, Colorado, 80305, United States
## 115 850, Hartford Drive, Boulder, Boulder County, Colorado, 80305, United States
## 200 759, Nighthawk Circle, Louisville, Boulder County, Colorado, 80027, United States
## zip_code city_update
## 1 80304 Boulder
## 2 80503 Longmont
## 3 80305 Boulder
## 4 80503 Niwot
## 5 80305 Boulder
## 6 80504 Longmont
## 7 80305 Boulder
## 8 80303 Boulder
## 10 80301 Boulder
## 11 80301 Boulder
## 13 80544 Niwot
## 14 80032 States
## 16 80305 Boulder
## 17 80306 Boulder
## 18 80302 Boulder
## 19 80021 States
## 20 80304 Boulder
## 21 80305 Boulder
## 22 80027 Louisville
## 23 80020 States
## 24 80302 Boulder
## 25 80027 Louisville
## 26 80304 Lane
## 28 80304 Boulder
## 29 80305 Boulder
## 31 80303 Boulder
## 34 80027 Boulder
## 35 80304 Boulder
## 36 80027 Louisville
## 37 80305 Boulder
## 38 80304 Boulder
## 39 80305 Boulder
## 40 80504 Longmont
## 41 80304 Lane
## 42 80301 Boulder
## 43 80303 Boulder
## 44 80025 Springs
## 45 80027 Superior
## 46 80304 Hills
## 47 80027 Superior
## 48 80027 Boulder
## 49 80305 Boulder
## 50 80304 Boulder
## 51 80305 Boulder
## 52 80301 Boulder
## 53 80027 Louisville
## 55 80305 Boulder
## 56 2872 Lafayette
## 57 80027 Louisville
## 58 80027 Louisville
## 59 80304 Boulder
## 60 80303 Boulder
## 61 80301 Boulder
## 62 80302 Boulder
## 63 80020 States
## 64 80301 Boulder
## 65 80304 Hill
## 66 80026 Erie
## 67 80026 Lafayette
## 69 80027 Louisville
## 70 80302 Boulder
## 71 80301 Boulder
## 72 80305 Boulder
## 73 80455 Acres
## 74 80026 Lafayette
## 75 80305 Boulder
## 76 80305 Boulder
## 78 80503 Niwot
## 80 80301 Boulder
## 81 80305 Boulder
## 82 80503 Court
## 83 80503 Court
## 84 80030 States
## 85 80305 Boulder
## 86 80503 Longmont
## 87 80304 Boulder
## 88 80304 Boulder
## 93 80021 States
## 94 80021 States
## 95 80023 States
## 97 80026 Lafayette
## 99 80305 Boulder
## 101 80303 Boulder
## 103 80301 Boulder
## 106 80305 Boulder
## 115 80305 Boulder
## 200 80027 Louisville
## loc_status
## 1 outside - Complete data throughout fire period
## 2 outside - Complete data throughout fire period
## 3 outside - Complete data throughout fire period
## 4 outside - Complete data throughout fire period
## 5 outside - Sensor offline during fire, returned online
## 6 outside - Sensor offline during fire, did not return online
## 7 outside - Sensor offline during fire, returned online
## 8 outside - Sensor offline during fire, returned online
## 10 outside - Complete data throughout fire period
## 11 inside - Complete data throughout fire period
## 13 inside - Complete data throughout fire period
## 14 outside - Complete data throughout fire period
## 16 outside - Complete data throughout fire period
## 17 outside - Complete data throughout fire period
## 18 outside - Complete data throughout fire period
## 19 outside - Complete data throughout fire period
## 20 inside - Complete data throughout fire period
## 21 outside - Sensor offline during fire, returned online
## 22 outside - Sensor offline during fire, returned online
## 23 outside - Complete data throughout fire period
## 24 outside - Complete data throughout fire period
## 25 outside - Sensor offline during fire, returned online
## 26 outside - Sensor offline during fire, returned online
## 28 outside - Complete data throughout fire period
## 29 inside - Complete data throughout fire period
## 31 outside - Complete data throughout fire period
## 34 inside - Sensor offline during fire, returned online
## 35 inside - Sensor offline during fire, returned online
## 36 outside - Complete data throughout fire period
## 37 inside - Complete data throughout fire period
## 38 outside - Complete data throughout fire period
## 39 outside - Complete data throughout fire period
## 40 outside - Complete data throughout fire period
## 41 inside - Sensor offline during fire, returned online
## 42 inside - Complete data throughout fire period
## 43 inside - Complete data throughout fire period
## 44 outside - Sensor offline during fire, returned online
## 45 outside - Sensor offline during fire, did not return online
## 46 outside - Complete data throughout fire period
## 47 outside - Sensor offline during fire, returned online
## 48 outside - Sensor offline during fire, returned online
## 49 outside - Complete data throughout fire period
## 50 outside - Complete data throughout fire period
## 51 inside - Complete data throughout fire period
## 52 inside - Complete data throughout fire period
## 53 inside - Complete data throughout fire period
## 55 outside - Sensor offline during fire, returned online
## 56 outside - Complete data throughout fire period
## 57 outside - Complete data throughout fire period
## 58 outside - Complete data throughout fire period
## 59 outside - Complete data throughout fire period
## 60 outside - Complete data throughout fire period
## 61 outside - Complete data throughout fire period
## 62 outside - Complete data throughout fire period
## 63 outside - Complete data throughout fire period
## 64 inside - Complete data throughout fire period
## 65 inside - Sensor offline during fire, returned online
## 66 outside - Complete data throughout fire period
## 67 outside - Complete data throughout fire period
## 69 inside - Sensor offline during fire, returned online
## 70 inside - Sensor offline during fire, did not return online
## 71 outside - Complete data throughout fire period
## 72 inside - Sensor offline during fire, returned online
## 73 inside - Sensor offline during fire, returned online
## 74 inside - Complete data throughout fire period
## 75 inside - Complete data throughout fire period
## 76 inside - Sensor offline during fire, returned online
## 78 outside - Complete data throughout fire period
## 80 outside - Complete data throughout fire period
## 81 inside - Sensor offline during fire, returned online
## 82 outside - Sensor offline during fire, returned online
## 83 outside - Sensor offline during fire, returned online
## 84 outside - Complete data throughout fire period
## 85 outside - Complete data throughout fire period
## 86 outside - Complete data throughout fire period
## 87 inside - Complete data throughout fire period
## 88 outside - Complete data throughout fire period
## 93 outside - Complete data throughout fire period
## 94 inside - Complete data throughout fire period
## 95 outside - Complete data throughout fire period
## 97 outside - Complete data throughout fire period
## 99 outside - Sensor offline during fire, returned online
## 101 outside - Complete data throughout fire period
## 103 outside - Complete data throughout fire period
## 106 inside - Sensor offline during fire, returned online
## 115 inside - Sensor offline during fire, returned online
## 200 inside - Sensor offline during fire, did not return online
# map of damaged/destroyed buildings and AQ sensors (indoor/outdoor)
# chose not to do icons for destroyed/damaged buildings because the map got really busy
(fire_destruction_and_AQ_plot = leaflet(sensors) %>%
addTiles() %>%
addAwesomeMarkers(icon = ~sensorIcon[Location],
lng = sensors$Lon, lat = sensors$Lat) %>%
addCircleMarkers(color = ~pal(destroyed_damaged$type),
radius = 3.5,
opacity = 1,
lng = destroyed_damaged$long, lat = destroyed_damaged$lat) %>%
addLegend("topright", pal = pal, values = ~destroyed_damaged$type,
title = "Fire damage"))
# # map based on when sensor was added
# (fire_destruction_and_AQ_plot = leaflet(sensors) %>%
# addTiles() %>%
# addAwesomeMarkers(icon = ~timeSensorIcon[loc_status],
# lng = sensors$Lon, lat = sensors$Lat) %>%
# addCircleMarkers(color = ~pal(destroyed_damaged$type),
# radius = 3.5,
# opacity = 1,
# lng = destroyed_damaged$long, lat = destroyed_damaged$lat) %>%
# addLegend("topright", pal = pal, values = ~destroyed_damaged$type,
# title = "Fire damage"))
# create pallete for location + time period
timeSensorIcon = awesomeIconList(
"inside - Complete data throughout fire period" = makeAwesomeIcon(
icon = "home",
iconColor = "white",
markerColor = "green",
library = "fa"
),
"inside - Sensor offline during fire, did not return online" = makeAwesomeIcon(
icon = "tree",
iconColor = "white",
markerColor = "gray",
library = "fa"
),
"inside - Sensor offline during fire, returned online" = makeAwesomeIcon(
icon = "home",
iconColor = "white",
markerColor = "blue",
library = "fa"
),
"outside - Complete data throughout fire period" = makeAwesomeIcon(
icon = "home",
iconColor = "white",
markerColor = "green",
library = "fa"
),
"outside - Sensor offline during fire, did not return online" = makeAwesomeIcon(
icon = "tree",
iconColor = "white",
markerColor = "gray",
library = "fa"
),
"outside - Sensor offline during fire, returned online" = makeAwesomeIcon(
icon = "tree",
iconColor = "white",
markerColor = "blue",
library = "fa"
))
# map based on when sensor was added
(fire_destruction_and_AQ_plot = leaflet(sensors) %>%
addTiles() %>%
addAwesomeMarkers(icon = ~timeSensorIcon[loc_status],
lng = sensors$Lon, lat = sensors$Lat) %>%
addCircleMarkers(color = ~pal(destroyed_damaged$type),
radius = 3.5,
opacity = 1,
lng = destroyed_damaged$long, lat = destroyed_damaged$lat) %>%
addLegend("topright", pal = pal, values = ~destroyed_damaged$type,
title = "Fire damage") %>%
addLegend("bottomright",
pal = colorFactor(c("green","gray","blue"), domain = unique(sensors$Status)),
values = ~sensors$Status, title = "Sensor Time Info"))
## Create spacetime objects -- 10 minute
# first, re-format datetime as a xts object
AQ_df$datetime = as.POSIXct(AQ_df$datetime)
# create spatial points object for each sensor
PA_sensors = SpatialPoints(AQ_df[!duplicated(AQ_df$ID), c("Lon", "Lat")],proj4string = CRS(prg))
summary(PA_sensors)
## Object of class SpatialPoints
## Coordinates:
## min max
## Lon -105.35808 -105.02234
## Lat 39.82646 40.13818
## Is projected: FALSE
## proj4string : [+proj=longlat +datum=WGS84 +no_defs]
## Number of points: 144
# contruct spatiotemporal object
PA_STFDF = stConstruct(AQ_df, space = c("Lon","Lat"), time = "datetime", crs = CRS(prg), SpatialObj = sensors)
#hrly_STFDF = stConstruct(AQ_df_hr, space = c("Lon","Lat"), time = "datetime", crs = CRS(prg), SpatialObj = stations)
# turn ST object into a STFDF (stationary points)
PA_STFDF = as(PA_STFDF,"STFDF")
#hrly_STFDF = as(hrly_STFDF,"STFDF")
# see summary of data
summary(PA_STFDF)
## Object of class STFDF
## with Dimensions (s, t, attr): (144, 17706, 12)
## [[Spatial:]]
## Object of class SpatialPoints
## Coordinates:
## min max
## Lon -105.35808 -105.02234
## Lat 39.82646 40.13818
## Is projected: FALSE
## proj4string : [+proj=longlat +datum=WGS84 +no_defs]
## Number of points: 144
## [[Temporal:]]
## Index timeIndex
## Min. :2021-12-30 00:00:00 Min. : 1
## 1st Qu.:2022-01-29 17:42:30 1st Qu.: 4427
## Median :2022-03-01 11:25:00 Median : 8854
## Mean :2022-03-01 11:25:00 Mean : 8854
## 3rd Qu.:2022-04-01 06:07:30 3rd Qu.:13280
## Max. :2022-05-01 23:50:00 Max. :17706
## [[Data attributes:]]
## X ID Name Location
## Min. : 1 Min. : 3062 Length:2549664 Length:2549664
## 1st Qu.: 482877 1st Qu.: 82879 Class :character Class :character
## Median : 965752 Median :110042 Mode :character Mode :character
## Mean : 965752 Mean : 99993
## 3rd Qu.:1448628 3rd Qu.:128801
## Max. :1931503 Max. :146952
## NA's :618161 NA's :618161
## pm25_a pm25_b temp rh
## Min. : 0.0 Min. : 0.0 Min. :-234.0 Min. : 0.0
## 1st Qu.: 0.0 1st Qu.: 0.0 1st Qu.: 38.2 1st Qu.:19.0
## Median : 0.5 Median : 0.4 Median : 50.0 Median :27.4
## Mean : 34.6 Mean : 6.6 Mean : 50.2 Mean :30.0
## 3rd Qu.: 3.9 3rd Qu.: 3.7 3rd Qu.: 62.6 3rd Qu.:40.4
## Max. :4999.0 Max. :3542.5 Max. : 127.4 Max. :94.6
## NA's :618161 NA's :618161 NA's :652435 NA's :652435
## a_b_agree pm_avg ab_difference corrected_pm
## Length:2549664 Min. : 0.0 Min. :-3537.2 Min. : -2.0
## Class :character 1st Qu.: 0.0 1st Qu.: -0.1 1st Qu.: 3.3
## Mode :character Median : 0.6 Median : 0.0 Median : 4.2
## Mean : 20.6 Mean : 28.0 Mean : 18.6
## 3rd Qu.: 4.0 3rd Qu.: 0.2 3rd Qu.: 5.4
## Max. :2514.5 Max. : 4994.6 Max. :1978.1
## NA's :618161 NA's :618161 NA's :652435
#summary(hrly_STFDF)
# object class
class(PA_STFDF)
## [1] "STFDF"
## attr(,"package")
## [1] "spacetime"
#class(hrly_STFDF)
# object dimensions
dim(PA_STFDF)
## space time variables
## 144 17706 12
#dim(hrly_STFDF)
ggplot(all_bounds, aes(geometry=geometry)) +
geom_sf() +
geom_sf(data=marshall_fire$geometry, fill="red", alpha=0.6) +
geom_sf(data=st_as_sf(PA_sensors))
## Time series plots
# non_null_stations = PA_STFDF[na.omit(stations)]
# temp and rh ts
stplot(PA_STFDF[,"2021-12-30::2021-12-31","temp"],mode="ts")
stplot(PA_STFDF[,"2021-12-30::2021-12-31","rh"],mode="ts")
stplot(PA_STFDF[,"2021-12-30::2021-12-31","corrected_pm"],mode="ts")
stplot(PA_STFDF[,"2021-12-30::2021-12-31","temp"],mode="xt")
stplot(PA_STFDF[,"2021-12-30::2021-12-31","rh"],mode="xt")
stplot(PA_STFDF[,"2021-12-30::2021-12-31","corrected_pm"],mode="xt")
# make an STFDF with only the sensors that had complete data for the fire
complete_sensors = sensors %>% filter(Status == "Complete data throughout fire period")
complete_STFDF = PA_STFDF[PA_STFDF@data$ID %in% complete_sensors$ID]
stplot(complete_STFDF[,,"corrected_pm"], mode="ts")
stplot(complete_STFDF[,"2021-12-30::2021-12-31","corrected_pm"], mode="ts")
stplot(complete_STFDF[,"2021-12-30::2021-12-31","corrected_pm"], mode="tp")
# Index 15 -> ID 60517, STANDLEY LAKE
plot(fortify(PA_STFDF[PA_STFDF@data$ID == 60517]["2021-12-30::2022-01-02", "corrected_pm"]))
# looks normal
# Index 23 -> ID 91877, Eisenhower Dr
plot(fortify(PA_STFDF[PA_STFDF@data$ID == 91877]["2021-12-30::2022-01-02", "corrected_pm"]))
# has incredibly large values midday Friday - Sunday, could be snow related?
# Index 36 -> ID 112606, 1333 Wildwood Ct., Boulder, CO 80305 sensor
plot(fortify(PA_STFDF[PA_STFDF@data$ID == 112606]["2021-12-30::2022-01-02", "corrected_pm"]))
# looking at purpleair map, A channel is 'downgraded', 0% confidence
# Index 40 -> ID 119547, Broadlands sensor
plot(fortify(PA_STFDF[PA_STFDF@data$ID == 119547]["2021-12-30::2022-01-02", "corrected_pm"]))
# reasonable values
Starting with Dec 30 - Jan 2, only outside sensors.
# reproject the data to utm?
reproj <- "+proj=utm +zone=13 +datum=WGS84 +units=km"
krig_data <- complete_STFDF[,"2021-12-30::2022-01-02"]
krig_data <- krig_data[krig_data@data$Location == "outside"]
# drop sensor 112606 b/c 0% confidence according to PurpleAir
krig_data <- krig_data[krig_data@data$ID != 112606]
# drop values for "Eisenhower Dr" sensor that're over 1000
krig_data@data[(krig_data@data$ID == 91877) &
(as.numeric(krig_data@data$corrected_pm) > 1000) &
!is.na(krig_data@data$corrected_pm),]$corrected_pm <- NA
# get rid of the weird sensor that might be throwing residuals off
krig_data <- krig_data[krig_data@data$Name != "ntsky"]
krig_data@sp <- spTransform(krig_data@sp, reproj)
krig_data <- aggregate(krig_data, by="hour", mean)
# split the data into the days we want
dec30 <- krig_data[, "2021-12-30"]
dec31 <- krig_data[, "2021-12-31"]
jan1 <- krig_data[, "2022-01-01"]
jan2 <- krig_data[, "2022-01-02"]
# split each day into the time periods of focus
# pre-fire period: betweeen 00:00 - 11:00 MST on 12/30/2021
pre_fire <- aggregate(dec30, by="11 hours", mean)[, "2021-12-30 00:00:00"]
# fire period (11am - 5pm)
during_fire <- aggregate(dec30[, "2021-12-30 11:00:00::2021-12-30 16:00:00"], by="6 hours", mean)
# post-fire on 12/30/2021
evening_of_fire <- aggregate(dec30[, "2021-12-30 16:00:00::2021-12-31 00:00:00"], by="8 hours", mean)
# dec 31st
december31 <- aggregate(dec31, by="day", mean)
# january 1st
january1 <- aggregate(jan1, by="day", mean)
# january 2nd
january2 <- aggregate(jan2, by="day", mean)
periods <- c(pre_fire, during_fire, evening_of_fire, december31, january1, january2)
# create variograms for these time periods
pre_fire.vgm <- variogram(corrected_pm~1, pre_fire[!is.na(pre_fire$corrected_pm),])
during_fire.vgm <- variogram(corrected_pm~1, during_fire[!is.na(during_fire$corrected_pm),])
evening_of_fire.vgm <- variogram(corrected_pm~1, evening_of_fire[!is.na(evening_of_fire$corrected_pm),])
december31.vgm <- variogram(corrected_pm~1, december31[!is.na(december31$corrected_pm),])
january1.vgm <- variogram(corrected_pm~1, january1[!is.na(january1$corrected_pm),])
january2.vgm <- variogram(corrected_pm~1, january2[!is.na(january2$corrected_pm),])
# plot the variograms
plot(pre_fire.vgm)
plot(during_fire.vgm)
plot(evening_of_fire.vgm)
plot(december31.vgm)
plot(january1.vgm)
plot(january2.vgm)
for (period in periods) {
cloud <- variogram(corrected_pm~1, period[!is.na(period$corrected_pm), ], cloud=TRUE)
map <- variogram(corrected_pm~1, period[!is.na(period$corrected_pm), ], map=TRUE, cutoff=20, width=2)
print(plot(cloud))
print(plot(map))
}
for (i in 0:11) {
if (i < 10) {
i <- paste("0", i, sep="")
}
t <- paste("2021-12-30 ", i, ":00:00", sep="")
i.data <- dec30[, t]
i.vgm <- variogram(corrected_pm~1, i.data[!is.na(i.data$corrected_pm),])
# print(plot(i.vgm, main=paste("Variogram at", t)))
i.fitted <- fit.variogram(i.vgm, vgm(c("Exp", "Sph", "Gau", "Mat")))
print(i.fitted)
print(plot(i.vgm, model=i.fitted, main=paste("Variogram at", t)))
}
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.000000 0.000000
## 2 Sph 3.904079 5.151113
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.000000 0.000000
## 2 Sph 3.095985 3.763504
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.000000 0.000000
## 2 Gau 3.260366 1.309588
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.4409823 0.000000
## 2 Sph 1.0171682 3.240835
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.1678131 0.00000
## 2 Gau 0.2917075 10.77585
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.0501557 0.00000
## 2 Gau 1.7184797 71.28375
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.028565 0.00000
## 2 Gau 0.201726 17.88329
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : linear model has singular covariance matrix
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : linear model has singular covariance matrix
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.01918062 0.00000
## 2 Sph 0.01776751 1.34037
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 0.034425382 0.00000
## 2 Sph 0.005219032 29.42497
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## model psill range
## 1 Nug 6181.04 0.0000
## 2 Exp 50155.36 -305.7408
## model psill range
## 1 Nug 6181.04 0.0000
## 2 Exp 50155.36 -305.7408
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## model psill range
## 1 Nug 391.9872 0.0000
## 2 Exp 375.1265 -189.0302
## model psill range
## 1 Nug 391.9872 0.0000
## 2 Exp 375.1265 -189.0302
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## model psill range
## 1 Nug 1.259568 0.000000
## 2 Gau 0.000000 4.613097
No good models variogram models!!!!
pre_fire.fitted <- fit.variogram(pre_fire.vgm, vgm(c("Exp", "Sph", "Gau", "Mat")))
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : singular model in variogram fit
pre_fire.fitted
## model psill range
## 1 Nug 86.24564 0.0
## 2 Exp 62387.14334 -10613.7
plot(pre_fire.vgm, model=pre_fire.fitted)
## model psill range
## 1 Nug 86.24564 0.0
## 2 Exp 62387.14334 -10613.7
during_fire.fitted <- fit.variogram(during_fire.vgm, vgm(c("Exp", "Sph", "Gau", "Mat")))
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
during_fire.fitted
## model psill range
## 1 Nug 152.5208 0.00000
## 2 Exp 938.8172 68.95887
plot(during_fire.vgm, model=during_fire.fitted)
evening_of_fire.fitted <- fit.variogram(evening_of_fire.vgm, vgm(c("Exp", "Sph", "Gau", "Mat")))
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, model, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
evening_of_fire.fitted
## model psill range
## 1 Nug 0.0000 0.00000
## 2 Gau 228.5243 23.99524
plot(evening_of_fire.vgm, model=evening_of_fire.fitted)
december31.fitted <- fit.variogram(december31.vgm, vgm(c("Exp", "Sph", "Gau", "Mat")))
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
december31.fitted
## model psill range
## 1 Nug 1.141338 0.00000
## 2 Gau 38.915184 16.75872
plot(december31.vgm, model=december31.fitted)
january1.fitted <- fit.variogram(january1.vgm, vgm(c("Exp", "Sph", "Gau", "Mat")))
## Warning in fit.variogram(object, x, fit.sills = fit.sills, fit.ranges =
## fit.ranges, : No convergence after 200 iterations: try different initial values?
january1.fitted
## model psill range
## 1 Nug 4.897684 0.000000
## 2 Sph 7.572974 2.524618
plot(january1.vgm, model=january1.fitted)
january2.fitted <- fit.variogram(january2.vgm, vgm(c("Exp", "Sph", "Gau", "Mat")))
january2.fitted
## model psill range
## 1 Nug 1.572919 0.00000
## 2 Exp 33.290787 12.57439
plot(january1.vgm, model=january2.fitted)
grd <- SpatialPixels(SpatialPoints(as_Spatial(st_make_grid(all_bounds, n=50))), proj4string = prg)
plot(grd)
grd <- grd[as_Spatial(all_bounds),]
plot(grd)
grd <- spTransform(grd, CRS(reproj))
plot(pre_fire$corrected_pm)
plot(during_fire$corrected_pm)
plot(evening_of_fire$corrected_pm)
plot(december31$corrected_pm)
plot(january1$corrected_pm)
plot(january2$corrected_pm)
pre_fire.idw <- idw(corrected_pm~1, pre_fire[!is.na(pre_fire$corrected_pm),], grd)
## [inverse distance weighted interpolation]
spplot(pre_fire.idw[,1], main="Pre-Fire Period on 12/30")
during_fire.idw <- idw(corrected_pm~1, during_fire[!is.na(during_fire$corrected_pm),], grd)
## [inverse distance weighted interpolation]
spplot(during_fire.idw[,1], main="During Fire on 12/30")
evening_of_fire.idw <- idw(corrected_pm~1, evening_of_fire[!is.na(evening_of_fire$corrected_pm),], grd)
## [inverse distance weighted interpolation]
spplot(evening_of_fire.idw[,1], main="After Fire on 12/30")
december31.idw <- idw(corrected_pm~1, december31[!is.na(december31$corrected_pm),], grd)
## [inverse distance weighted interpolation]
spplot(december31.idw[,1], main="December 31st")
january1.idw <- idw(corrected_pm~1, january1[!is.na(january1$corrected_pm),], grd)
## [inverse distance weighted interpolation]
spplot(january1.idw[,1], main="January 1st")
january2.idw <- idw(corrected_pm~1, january2[!is.na(january2$corrected_pm),], grd)
## [inverse distance weighted interpolation]
spplot(january2.idw[,1], main="January 2nd")
aqi_colors <- c("#00e400",
"#ffff00",
"#ff7e00",
"#ff0000",
"#8f3f97",
"#7e0023")
aqi_labels <- c("Good (0-12.0)",
"Moderate (12.1-35.4)",
"Unhealthy for Sensitive Groups (35.5-55.4)",
"Unhealthy (55.5-150.4)",
"Very Unhealthy (150.5-250.4)",
"Hazardous (250.5-500.4)")
factorize_aqi <- function (idw) {
idw$PM2.5 <- cut(idw$var1.pred, breaks=c(0, 12.0, 35.4, 55.4, 150.4, 250.4), na.omit=T)
return(idw)
}
plot_idw <- function(idw, title) {
return(ggplot(all_bounds) +
geom_sf() +
geom_sf(data=st_as_sf(factorize_aqi(idw)), aes(col=PM2.5), size=2) +
scale_color_manual(values=aqi_colors, labels=aqi_labels, name="PM2.5 (µg/m^3)") +
geom_sf(data=marshall_fire$geometry, col="gray20", alpha=0.1) +
ggtitle(title))
}
plot_idw(pre_fire.idw, "Pre Fire IDW")
plot_idw(during_fire.idw, "During Fire IDW")
plot_idw(evening_of_fire.idw, "Evening of fire IDW")
plot_idw(december31.idw, "December 31st IDW")
plot_idw(january1.idw, "January 1st IDW")
plot_idw(january2.idw, "January 2nd IDW")
pre_fire.idw_cv <- krige.cv(corrected_pm~1, pre_fire[!is.na(pre_fire$corrected_pm),], nfold=38)
pre_fire.idw_sd <- sd(pre_fire.idw_cv$residual)
paste("Pre fire standard deviation:", pre_fire.idw_sd)
## [1] "Pre fire standard deviation: 7.58649993004707"
spplot(pre_fire.idw_cv, main="Pre-Fire Period on 12/30")
pre_fire.idw_cv$zscore <- pre_fire.idw_cv$residual / pre_fire.idw_sd
during_fire.idw_cv <- krige.cv(corrected_pm~1, during_fire[!is.na(during_fire$corrected_pm),], nfold=38)
during_fire.idw_sd <- sd(during_fire.idw_cv$residual)
paste("During fire standard deviation:", during_fire.idw_sd)
## [1] "During fire standard deviation: 15.376531044756"
spplot(during_fire.idw_cv, main="During Fire on 12/30")
during_fire.idw_cv$zscore <- during_fire.idw_cv$residual / during_fire.idw_sd
evening_of_fire.idw_cv <- krige.cv(corrected_pm~1, evening_of_fire[!is.na(evening_of_fire$corrected_pm),], nfold=38)
evening_of_fire.idw_sd <- sd(evening_of_fire.idw_cv$residual)
paste("Evening of fire standard deviation:", evening_of_fire.idw_sd)
## [1] "Evening of fire standard deviation: 4.85088863186492"
spplot(evening_of_fire.idw_cv, main="After Fire on 12/30")
evening_of_fire.idw_cv$zscore <- evening_of_fire.idw_cv$residual / evening_of_fire.idw_sd
december31.idw_cv <- krige.cv(corrected_pm~1, december31[!is.na(december31$corrected_pm),], nfold=38)
december31.idw_sd <- sd(december31.idw_cv$residual)
paste("December 31st standard deviation:", december31.idw_sd)
## [1] "December 31st standard deviation: 2.36672112753216"
spplot(december31.idw_cv, main="December 31st")
december31.idw_cv$zscore <- december31.idw_cv$residual / december31.idw_sd
january1.idw_cv <- krige.cv(corrected_pm~1, january1[!is.na(january1$corrected_pm),], nfold=38)
january1.idw_sd <- sd(january1.idw_cv$residual)
paste("January 1 standard deviation:", january1.idw_sd)
## [1] "January 1 standard deviation: 3.68472999489718"
spplot(january1.idw_cv, main="January 1st")
january1.idw_cv$zscore <- january1.idw_cv$residual / january1.idw_sd
january2.idw_cv <- krige.cv(corrected_pm~1, january2[!is.na(january2$corrected_pm),], nfold=38)
january2.idw_sd <- sd(january2.idw_cv$residual)
paste("January 2 standard deviation:", january2.idw_sd)
## [1] "January 2 standard deviation: 3.6116917414784"
spplot(january1.idw_cv, main="January 2nd")
january2.idw_cv$zscore <- january2.idw_cv$residual / january2.idw_sd
plot_idw_resid <- function(idw.cv, title) {
return(ggplot(all_bounds) +
geom_sf(fill="#FbFbFb", color="gray") +
geom_sf(data=marshall_fire$geometry, col="red", alpha=0.1) +
geom_sf(data=st_as_sf(idw.cv), aes(fill=residual), color="black", shape=21, size=2) +
scale_color_gradient2(aesthetics="fill") +
ggtitle(title))
}
plot_idw_resid_zscore <- function(idw.cv, title) {
return(ggplot(all_bounds) +
geom_sf(fill="#FbFbFb", color="gray") +
geom_sf(data=marshall_fire$geometry, col="red", alpha=0.1) +
geom_sf(data=st_as_sf(idw.cv), aes(fill=zscore), color="black", shape=21, size=2) +
scale_color_gradient2(aesthetics="fill") +
ggtitle(title))
}
plot_idw_cv <- function(idw, idw_title, idw.cv) {
p1 <- plot_idw(idw, idw_title)
p2 <- plot_idw_resid(idw.cv, "Residuals")
return(grid.arrange(p1, p2, ncol=2))
}
plot_idw_resid(pre_fire.idw_cv, "Pre-Fire Residuals")
plot_idw_resid(during_fire.idw_cv, "During Fire Residuals")
plot_idw_resid(evening_of_fire.idw_cv, "After Fire Residuals")
plot_idw_resid(december31.idw_cv, "December 31st Residuals")
plot_idw_resid(january1.idw_cv, "January 1st Residuals")
plot_idw_resid(january2.idw_cv, "January 2nd Residuals")
plot_idw_resid_zscore(pre_fire.idw_cv, "Pre-Fire Residual Z-Scores")
plot_idw_resid_zscore(during_fire.idw_cv, "During Fire Residual Z-Scores")
plot_idw_resid_zscore(evening_of_fire.idw_cv, "After Fire Residual Z-Scores")
plot_idw_resid_zscore(december31.idw_cv, "December 31st Residual Z-Scores")
plot_idw_resid_zscore(january1.idw_cv, "January 1st Residual Z-Scores")
plot_idw_resid_zscore(january2.idw_cv, "January 2nd Residual Z-Scores")
plot_idw_cv(pre_fire.idw, "Pre-Fire IDW", pre_fire.idw_cv)
plot_idw_cv(during_fire.idw, "During Fire IDW", during_fire.idw_cv)
plot_idw_cv(evening_of_fire.idw, "Evening of Fire IDW", evening_of_fire.idw_cv)
plot_idw_cv(december31.idw, "December 31st IDW", december31.idw_cv)
plot_idw_cv(january1.idw, "January 1st IDW", january1.idw_cv)
plot_idw_cv(january2.idw, "January 2nd IDW", january2.idw_cv)
# One sensor is "Eisenhower Dr"
plot(fortify(PA_STFDF[PA_STFDF@data$Name == "Eisenhower Dr"]["2021-12-30::2022-01-02", "corrected_pm"]))
# The other is "Purple House"
plot(fortify(PA_STFDF[PA_STFDF@data$Name == "Purple House"]["2021-12-30::2022-01-02", "corrected_pm"]), col="red")
The issue is coming from the Eisenhower Dr sensor most likely, as this is the one that seems messed up by the snow. Going to see how everything looks excluding the high values for this sensor.
plot(fortify(PA_STFDF[PA_STFDF@data$Name == "Eisenhower Dr"]["2021-12-30::2022-01-02", "corrected_pm"]) %>%
filter(corrected_pm <= 1000))
Checking the sensor on the edge of the fire boundary with high residuals/z scores.
plot(fortify(PA_STFDF[PA_STFDF@data$Name == "ntsky"]["2021-12-30::2022-01-02", "corrected_pm"]))
for (j in c(30, 31, 01, 02)) {
if (j > 3) {
date <- paste("2021-12-", j, sep="")
} else {
date <- paste("2022-01-", j, sep="")
}
for (i in 0:23) {
if (i < 10) {
i <- paste("0", i, sep="")
}
t <- paste(date, " ", i, ":00:00", sep="")
i.data <- krig_data[, t]
i.idw <- idw(corrected_pm~1, i.data[!is.na(i.data$corrected_pm),], grd)
print(plot_idw(i.idw, paste(t, "IDW")))
}
}
for (j in c(30, 31, 01, 02)) {
if (j > 3) {
date <- paste("2021-12-", j, sep="")
} else {
date <- paste("2022-01-", j, sep="")
}
for (i in 0:23) {
if (i < 10) {
i <- paste("0", i, sep="")
}
t <- paste(date, " ", i, ":00:00", sep="")
i.data <- krig_data[, t]
i.idw <- idw(corrected_pm~1, i.data[!is.na(i.data$corrected_pm),], grd)
i.idw_cv <- krige.cv(corrected_pm~1, i.data[!is.na(i.data$corrected_pm),], nfold=38)
i.idw_sd <- sd(i.idw_cv$residual)
i.idw_cv$zscore <- i.idw_cv$residual / i.idw_sd
print(plot_idw_resid_zscore(i.idw_cv, paste(t, "IDW Residual Z-Scores")))
}
}
# get the full timerange of data & clean it up
daily <- complete_STFDF[,"2022-01-01::"]
daily <- daily[daily@data$Location == "outside"]
# drop sensor 112606 b/c 0% confidence according to PurpleAir
daily <- daily[daily@data$ID != 112606]
# drop values for "Eisenhower Dr" sensor that're over 1000
daily@data[(daily@data$ID == 91877) &
(as.numeric(daily@data$corrected_pm) > 1000) &
!is.na(daily@data$corrected_pm),]$corrected_pm <- NA
# get rid of the weird sensor that might be throwing residuals off
daily <- daily[daily@data$Name != "ntsky"]
# get rid of "ponderosa" sensor values after april 11 (when sensor looks to get wonky)
daily@data[(daily@data$ID == 102614) & (daily@endTime < "2022-04-11") & !is.na(daily@data)] <- NA
# reproject & aggregate
daily@sp <- spTransform(daily@sp, reproj)
daily <- aggregate(daily, by="day", mean)
for (j in 1:4) {
date <- paste("2022-0", j, sep="")
if (j %in% c(1, 3, 5)) {
range = 1:31
} else if (j == 2) {
range = 1:28
} else {
range = 1:30
}
for (i in range) {
if (i < 10) {
i <- paste("0", i, sep="")
}
day <- paste(date, "-", i, sep="")
data <- daily[, day]
idw <- idw(corrected_pm~1, data[!is.na(data$corrected_pm),], grd)
print(plot_idw(idw, paste(day, "IDW")))
}
}
plot(daily[, "2022-04-25"]$corrected_pm) # sensor at index 29 is weird = sensor ID 102614
plot(fortify(PA_STFDF[PA_STFDF@data$ID == 102614]['2022-04', "corrected_pm"]))
# ponderosa sensor has 0% confidence on purpleair map at present date
# OLD HELPER FUNCTION USED TO DETERMINE LOCATION
get_sensor_info <- function(id) {
# for each id, get the rest of the information (including lat/long)
info = id_key[id_key$ID == id,]
# create a spatial version of info
info_spatial = st_as_sf(info, coords = c("Lon","Lat"), crs = prg)
#intersect ID info with the county map to determine location
city_intersect = st_intersection(info_spatial,st_buffer(fire_counties,0))["ZONEDESC"]
city = city_intersect$ZONEDESC
info = info %>%
mutate(city = ifelse(length(city) == 0, NA, city))
return(city)
}
# create a key coordinating the sensor ID number to its lat/lon
(id_key = AQ_df %>%
group_by(ID) %>%
dplyr::select(ID, Lon, Lat, Location) %>%
unique() %>%
# remove NAs for spatial intersections to occur
na.omit())
## # A tibble: 144 x 4
## # Groups: ID [144]
## ID Lon Lat Location
## <int> <dbl> <dbl> <chr>
## 1 3062 -105. 40.0 outside
## 2 4583 -105. 40.1 outside
## 3 4643 -105. 40.0 outside
## 4 15663 -105. 40.1 outside
## 5 27120 -105. 40.0 outside
## 6 30743 -105. 40.1 outside
## 7 33109 -105. 40.0 outside
## 8 38623 -105. 40.0 outside
## 9 45325 -105. 40.1 outside
## 10 50133 -105. 40.1 inside
## # ... with 134 more rows
for (sensor in 1:30) {
id <- id_key$ID[sensor]
sensor_info <- get_sensor_info(id)
# check that the sensor has some data for our time range
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# print("1")
if(length(sensor_xts$ID) > 0) {
# plot pm2.5 on one graph for both channels
# sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
# sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
# print("2")
p <- ggplot(fortify(sensor_xts["2021-12-30::2021-12-31", "corrected_pm"]),
aes(x=Index, y=as.numeric(corrected_pm))) +
geom_point() +
# geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
# scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="Corrected PM2.5",
title = paste(sensor_info, "-", "sensor", id, "\nCorrected PM2.5"))
print(p)
# print("3")
# png(file=paste("images/sensor", id, "_pm25_ts.png"), width=500, height=500)
# print(p)
# dev.off()
# plot temperature and relative humidity on one plot !!! this doesn't seem possible with ggplot
# sensor_temp <- sensor_xts["2021-12-30::2021-12-31", "temp"]
# sensor_rh <- sensor_xts["2021-12-30::2021-12-31", "rh"]
#
# p2 <- ggplot(fortify(sensor_temp), aes(x=Index, y=temp)) +
# geom_line() +
# labs(x="date",
# title = paste(sensor_info["city"], "-", sensor_info["Location"],
# "\nsensor", id, "temperature"))
#
# p3 <- ggplot(fortify(sensor_rh), aes(x=Index, y=rh)) +
# geom_line() +
# labs(x="date",
# title = paste(sensor_info["city"], "-", sensor_info["Location"],
# "\nsensor", id, "relative humidity"))
#
# plots <- gridExtra::grid.arrange(p, arrangeGrob(p2, p3, ncol=2), nrow=2)
# print(plots)
# png(file=paste("images/sensor", id, "_ts.png"), width=500, height=500)
# print(plots)
# dev.off()
}
}
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
## dist is assumed to be in decimal degrees (arc_degrees).
## although coordinates are longitude/latitude, st_intersection assumes that they are planar
# save output to enhance resolution
jpeg(file = "../images/timeplot.jpeg")
stplot(PA_STFDF[,"2021-12-30::2021-12-31","corrected_pm"])
dev.off()
## png
## 2
stplot(PA_STFDF[,,"ID"], mode="ts")
Pairs gotten by looking at the purpleair sensor map, not every sensor has a close one, so I just picked the closest outdoor one we have data for. Generally on the same street, just a few houses down.
# sensor 103724 - Louisville inside
id <- 103724
sensor_info <- get_sensor_info(id)
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# plot pm2.5 on one graph for both channels
sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
p <- ggplot(fortify(sensor_pm25a), aes(x=Index, y=pm25_a, color="black")) +
geom_line() +
geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="PM2.5",
title = paste(sensor_info$city, "-", sensor_info$Location, sensor_info$Lon,",", sensor_info$Lat,
"\nsensor", id, "pm2.5"))
# sensor 134408 - Louisville outside
id <- 65669
sensor_info <- get_sensor_info(id)
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# plot pm2.5 on one graph for both channels
sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
p1 <- ggplot(fortify(sensor_pm25a), aes(x=Index, y=pm25_a, color="black")) +
geom_line() +
geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="PM2.5",
title = paste(sensor_info$city, "-", sensor_info$Location, sensor_info$Lon,",", sensor_info$Lat,
"\nsensor", id, "pm2.5"))
plots <- gridExtra::grid.arrange(p, p1, nrow=2)
print(plots)
# sensor 141956 - Louisville inside
id <- 141956
sensor_info <- get_sensor_info(id)
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# plot pm2.5 on one graph for both channels
sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
p <- ggplot(fortify(sensor_pm25a), aes(x=Index, y=pm25_a, color="black")) +
geom_line() +
geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="PM2.5",
title = paste(sensor_info$city, "-", sensor_info$Location, sensor_info$Lon,",", sensor_info$Lat,
"\nsensor", id, "pm2.5"))
# sensor 81149 - Louisville outside
id <- 81149
sensor_info <- get_sensor_info(id)
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# plot pm2.5 on one graph for both channels
sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
p1 <- ggplot(fortify(sensor_pm25a), aes(x=Index, y=pm25_a, color="black")) +
geom_line() +
geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="PM2.5",
title = paste(sensor_info$city, "-", sensor_info$Location, sensor_info$Lon,",", sensor_info$Lat,
"\nsensor", id, "pm2.5"))
plots <- gridExtra::grid.arrange(p, p1, nrow=2)
print(plots)
# sensor 79097 - Louisville inside
id <- 79097
sensor_info <- get_sensor_info(id)
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# plot pm2.5 on one graph for both channels
sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
p <- ggplot(fortify(sensor_pm25a), aes(x=Index, y=pm25_a, color="black")) +
geom_line() +
geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="PM2.5",
title = paste(sensor_info$city, "-", sensor_info$Location, sensor_info$Lon,",", sensor_info$Lat,
"\nsensor", id, "pm2.5"))
# sensor 81149 - Lousiville outside
id <- 81149
sensor_info <- get_sensor_info(id)
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# plot pm2.5 on one graph for both channels
sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
p1 <- ggplot(fortify(sensor_pm25a), aes(x=Index, y=pm25_a, color="black")) +
geom_line() +
geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="PM2.5",
title = paste(sensor_info$city, "-", sensor_info$Location, sensor_info$Lon,",", sensor_info$Lat,
"\nsensor", id, "pm2.5"))
plots <- gridExtra::grid.arrange(p, p1, nrow=2)
print(plots)
# sensor 88533 - Louisville inside
id <- 88533
sensor_info <- get_sensor_info(id)
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# plot pm2.5 on one graph for both channels
sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
p <- ggplot(fortify(sensor_pm25a), aes(x=Index, y=pm25_a, color="black")) +
geom_line() +
geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="PM2.5",
title = paste(sensor_info$city, "-", sensor_info$Location, sensor_info$Lon,",", sensor_info$Lat,
"\nsensor", id, "pm2.5"))
# sensor 92321 - Louisville outside
id <- 92321
sensor_info <- get_sensor_info(id)
sensor_xts <- PA_STFDF[PA_STFDF@data$ID == id]
# plot pm2.5 on one graph for both channels
sensor_pm25a <- sensor_xts["2021-12-30::2021-12-31", "pm25_a"]
sensor_pm25b <- sensor_xts["2021-12-30::2021-12-31", "pm25_b"]
p1 <- ggplot(fortify(sensor_pm25a), aes(x=Index, y=pm25_a, color="black")) +
geom_line() +
geom_line(data=fortify(sensor_pm25b), aes(y=pm25_b), color="tomato") +
scale_color_manual(name = "Sensor Data", values = c("PM2.5 A Channel" = "black", "PM2.5 B Channel" = "tomato")) +
labs(x="date",
y="PM2.5",
title = paste(sensor_info$city, "-", sensor_info$Location, sensor_info$Lon,",", sensor_info$Lat,
"\nsensor", id, "pm2.5"))
plots <- gridExtra::grid.arrange(p, p1, nrow=2)
print(plots)
# why is one showing up way off -- investigate, I wonder if this sensor moved during the process?
# sensor 112606
# png(file = "images/fire_map2.png", width = 4500, height = 3800, res = 300)
# (detailed_sensor_map = basemap +
# geom_point(data=(plot_data %>% filter(ID != "112606")), aes(x=Lon,y=Lat, label = ID, shape = Location, color = Status, label = rownames(ID))) +
# scale_colour_brewer(palette = "Set1") +
# theme_light() +
# theme(legend.position = "bottom", legend.box = "vertical") +
# #guides(scale = FALSE) +
# ggtitle("Map of sensors in Louisville and Superior, 12-30-21 to 2-26-22"))
# dev.off()
# add sensors, check which sensor is off
# (ggplot(plot_data, aes(x=Lon,y=Lat, label = ID, shape = Location, color = fire_period, label = ID))) +
# geom_point() +
# geom_text(check_overlap = TRUE)
# Create time period variable based on availability during fire, inside/outside status for mapping
# (plot_data = AQ_df %>%
# group_by(ID) %>%
# dplyr::select(ID, Lon, Lat, Location) %>%
# # use location[1] presuming indoor/outdoor doesn't change
# dplyr::summarize(Location = Location[1],
# Lat = Lat[1],
# Lon = Lon[1]))
# png(file = "images/month_map.png", width = 6000, height = 3800, res = 300)
# (month_plot = basemap +
# geom_point(data = month_added_plot, aes(x=Lon,y=Lat, label = ID, shape = Location, color = Month)) +
# facet_grid(cols = vars(Month)))
# dev.off()